High-Resolution Variograms from Optical Borehole Imaging for Borehole-Scale Geostatistical Rock Simulation and Lattice Boltzmann Flow Simulation in High-Permeability Carbonates; #120091 (2013)

Geometry of the pore space is a fundamental control on flow in carbonates over reservoir exploitation time scales. Fluid properties and pressure gradients are temporally variable and can be changed by use of effective reservoir production, workover, and stimulation procedures. The geometry of the porous medium largely remains constant except if significant pressure variations deform or fracture the reservoir or acid fracturing is applied. Formulating methods to accurately simulate subsurface carbonate pore systems is a difficult challenge. A technique is presented that uses spatial statistics applied to digital optical borehole images to create realizations of carbonate megaporosity. Megapores are “equant to equant-elongate pores whose average diameter is larger than 4 mm, and for tubular or platy pores whose average cross-sectional diameter or thickness, respectively, is larger than 4 mm” (Choquette and Pray, 1970). Flow within megapores at the borehole scale was simulated using Lattice Boltzmann methods. The borehole-scale Lattice Boltzmann computations are orders of magnitude larger than conventional pore-scale applications. Traditional geostatistical models are often based on sparse data obtained from cores or down-hole tests. Conversely, image data are dense and provide highly constrained variograms that reflect the detailed structure of the rock megapore system and its vertical arrangement within highfrequency depositional cycles. The high-resolution variograms obtained from borehole image data justify the fitting of unusually complex variogram models containing relevant stratigraphic information. Moreover, with and without conditioning to the borehole image data, satisfying stochastic realizations of the rock can be generated. Current research efforts focus on accurate quantification of permeability and non-Darcy flow parameter measurements in the karst Pleistocene carbonate rocks of the Biscayne aquifer in southeastern Florida, which can have exceptionally high porosity and permeability values. Standard laboratory core permeability measurement techniques (e.g., air permeameter) are not suited to accurately determine the permeability of these rocks. Furthermore, benchtop analyses of limestone core representative of extremely high megaporosity have never been accomplished due to breakage of relatively fragile vuggy-rock intervals during drilling. In the southeast Florida study area, extremely megaporous aquifer intervals are attributed to a megapore system formed by Ophiomorpha ichnofabric (Droser and Bottjer, 1989; Cunningham et al., 2009; Cunningham et al., 2010; Cunningham and Sukop, 2011). In one borehole representative of the Biscayne aquifer, imagery is available at approximately 2 mm/pixel resolution over 17 m of borehole depth. These data were combined with an interpolated caliper log to provide approximately 6 million available points on a variable-radius, 3-D representation of the borehole wall where each pixel was categorized as either pore space or solid for use in variogram computation. Sonic caliper data may enhance the ability of variograms to detect details of the spatial structure from the digital optical borehole images and enable higher fidelity rock simulations. The borehole image data would form approximately 10 unique pairs of points for variogram estimation producing a data set so large that variogram computation is not manageable with standard software (SGeMS; Remy et al., 2009; Bianchi and Zheng, 2009) on a conventional PC. Computations on 1, 5, and 20% subsamples of the 6-million-point data set were manageable and show that the variograms are well-resolved even at the 1% sampling rate. Three nested (superimposed) variogram models closely fit the horizontal data-based variogram with only slight anisotropy in the horizontal direction. The vertical variograms show strong short-range autocorrelation out to a lag of about 0.5 m and a number of superimposed periodicities (with 0.53and 3-m periods) to lags of approximately 8 m. These reflect megaporosity changes that occur within m-scale, intracycle high-frequency vertical lithofacies in the study area (Cunningham and Sukop, 2011). Six nested variogram models were fitted to the data-based vertical variogram estimates. A preliminary simulation of a rectangular 0.4 × 0.4 × 17 m rock prism surrounding the borehole has been completed (Figure 1). The rock simulation was 80 × 80 × 3400 voxels corresponding to a resolution of 5 mm per voxel. Sequential indicator simulation using the SISIM program of GSLIB (e.g., Remy et al., 2009; Deutch and Journel, 1998) conditioned to the observations was employed. It was necessary to modify the SISIM program and the subroutine COVA3 to: accommodate the large number of nested variogram models; simultaneously use hole-effect and damped hole-effect variogram models; and permit Gaussian models to be part of the variogram structure input to SISIM. As with variogram estimation, there are many possible permutations of simulation parameters and approaches that are potentially valid that can affect the results of the simulation. Moreover, the simulation process is stochastic but only a single realization has been considered due to the computational burden of the Lattice Boltzmann Methods (LBM) with limited processing power and lengthy run times. Ideally, an ensemble of realizations would be prepared and LBM flow simulation would be conducted in each of them. A decision to not pursue these computations was based strictly on computer limitations that are expected to continue to be alleviated as more computing power is routinely available. Use of the full dense borehole data for conditioning should appreciably reduce variability among rock realizations, and thereby limit the variability in the final flow simulation results. The rock prism was simulated as solid or pore at each voxel and the target proportion of pores and solids was specified. A zonal anisotropy approach (Deutch and Journel, 1998) was adopted to accommodate the different horizontal and vertical variograms. Lattice Boltzmann (LB) methods represent an alternative approach to simulate single and multiphase fluid flows. One reason for the popularity of LB methods among porous media flow modelers is the relative ease of incorporating the complex boundary condition of the pore/rock interface. Since Succi et al. (1989) first demonstrated LB estimation of the permeability in complex porous media, LB simulat ions have been widely used but at much smaller scales than those we consider here [e.g., total simulation domain size in parentheses, Ferréol and Rothman, 1995 (0.2 cm); Auzerais et al, 1996 (~0.2 cm); Pan et al, 2001 (~0.6 cm); Fredrich et al., 2006 (~0.3 cm); Wright et al., 2006 (~0.4 mm); Okabe and Blunt, 2007 (~0.7 mm); Polacci et al., 2009 (~0.2 cm); Boek and Venturolia, 2010 (0.25 cm); Degruyter et al., 2010 (~0.3 cm)]. For the purpose of estimating the permeability in the borehole region, the key result of LB simulation is the Darcy flux q that occurs through the rock prism under the applied pressure gradient. The Darcy flux can be obtained as the summation of the x-axis velocities over any plane perpendicular to the bulk flow direction divided by the bulk area. Darcy’s Law can be rearranged to solve for intrinsic permeability (k) when all other terms are known. In preliminary simulations, for example, the intrinsic permeability k of a simulated rock prism in LB length units (lu) was 0.21 lu. The following scaling relationship allows calculation of the corresponding physical intrinsic permeability: 2        LBM physical LBM physical L L k k . (1) The physical intrinsic permeability of a preliminary simulated rock prism in which the length scale ratio is 0.005 m/lu is 5.2 × 10 m or about 5,300,000 darcys. The preliminary Lattice Boltzmann-derived permeability for a rock prism agrees well with previous estimates computed on ~0.1-m cubic samples of the same aquifer materials using computed tomography scans (Cunningham et al., 2009; Cunningham et al., 2010, Cunningham and Sukop, 2011). The magnitude of the estimated permeability is also consistent with Poiseuille equation-based values of the permeability of smooth pipes of size comparable to the megapores in the simulated rock. For example, a smooth, straight 0.02 m diameter pipe has an intrinsic permeability of 1.25 × 10 m (12,700,000 darcys). Rock pores are not smooth pipes however, and roughness, tortuosity, imperfect connectivity, and the presence of low conductivity rock matrix will all reduce the permeability values. Additional refinement of the rock simulations and much larger scale simulations than the 17-m prisms considered here will be possible as computation capacities continue to grow. Furthermore, Lattice Boltzmann methods offer powerful multiphase-fluid simulation capabilities that can readily be applied to the simulated rocks.